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ABSTRACT 

In this work we employ an MHD numerical code to reproduce the morphology observed 
for SN 1006 in radio synchrotron and thermal X-ray emission. We introduce a density 
discontinuity, in the form of a flat cloud parallel to the Galactic Plane, in order to 
explain the NW filament observed in optical wavelengths and in thermal X-rays. We 
compare our models with observations. We also perform a test that contrasts the radio 
emitting bright limbs of the SNR against the central region, finding additional support 
to our results. Our main conclusion is that the most probable direction of the ambient 
magnetic field is on average perpendicular to the Galactic Plane. 

Key words: MHD - methods: numerical - ISM: supernova remnants - supernovae: 
individual: SN 1006 - radiation mechanism: general. 



1 INTRODUCTION 

SN 1006 is a young supernova remnant (SNR) whose mor- 
phology in radio-continuum can be described as a partial 
shell with a radius of 15^ Since most of the emission arises 
in two main bright arcs towards the NE and SW, this 
] SNR has been classified as bilateral or bar rel-shaped (BSNR 

■ iKesteven h Caswelll [TqStI : iGaenslerlllQQSh . Another charac- 

■ teristic of this source is that the thermal X-ray emission 
is uniformly distributed throughout the SNR, and is as- 

^ socia t ed to the ejected material (see Cassam-Chena'i et al.l 
I2OO8I : iMiceli et al.ll2009l ). Recentlv. IXcero et al.1 (|2007l ) car- 
ried out observations of this object with XMM-Newton, con- 
firming; the presen ce of a bright filament in the NW (see 
iDubner et al]l2002h . suggesting that the expanding shell is 
colliding into a denser medium. Non-thermal X-ray emission 
is predominant in the NE and SW almost coincident wit h 
the synchrotron radio em i ssion (|Cassam-Chena'i et al.ll2008h . 

iKesteven h Caswelll (|l987l ) propose two possible mech- 
anisms for the production of the radio-arcs: an intrinsic one, 
due for example to the interaction of internal jets with the 
SNR shell, and an extrinsic one, dominated by the charac- 
teristics of the surrounding ISM. 

Several theoretical works have lately studied the ex- 
trinsic mechanism for the formation of two opposed ra- 
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dio arcs in BSNRs through particle injection models to- 
g:ether with analytical sel f -similar (for recent works see 
ICassam-Chenai et al.1 12OO8I: IPetruk et al.ll2009l ) or magne- 
tohydrodynamic (see lOrlando et al.l 12007 ) numerical simu- 
lations. 



lOrlando et al.1 (l2007f ) explores, through MHD simula- 
tions, different scenarios for explaining the asymmetries be- 
tween the arcs in BSNRs. They propose two main mecha- 
nisms: the SNR expands in a non-uniform medium with an 
almost uniform magnetic field, and the SNR evolves in a uni- 
form ISM with a non-uniform magnetic field. In either case, 
the density gradient or the magnetic field gradient are not 
aligned with the line of sight. To generate the synthetic syn- 
chrotron maps, they take into account three possible ways 
of particle injection: isotropic, quasi-parallel, and quasi- 
perpendicular. The quasi-parallel (or perpendicular) nature 
of the injection is determined by the angle between the shock 
no rmal and the s hocke d magnetic field ((\)bs)- According 
to lOrlando et aD (l200/t ) only the quasi-perpendicular and 
isotropic injection models reproduce radio images similar to 
the o bserved ones (in agreement with iFulbright h Revnoldsl 
I1990I ). However it is important to know the direction of the 
unperturbed interstellar magnetic field (ISMF), in order to 
decide which injection mechanism is predominant. 



iRothenflug et al.1 (|2004l ) carried out an analysis of the 
radio-continuum and non-thermal X-ray emission of the 
remnant of SN 1006. By using an i?7r/3 criterion and based 
mainly on the non-th ermal X-ray emission distribution, 
IRothenflug et al.1 (|2004h conclude that the ISMF around SN 
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1006 is mainly perpendicular to the NE and SW bright radio 

arcs (i.e. the ISMF is parallel to the Galactic plane) . 

A different conclusion was obtained bv IPetruk et al.l 

(Eooi), who employed the analytical Sedov density and pres- 
sure profiles, where the aforementioned injection particles 
and self-similar magnetic field configuration are assumed to 
depend only on the shock obliquity through the compres- 
sion factor in order to model the radio-continuum emission 
from SNRs. They analyse the orientation of the ISMF based 
on the azimuthal profile of the radio surface brightness on 
both observed and simulated radio maps of SN 1006. Ap- 
plying this method, they found that the models that best 
reproduce both the observed radio morphology and the az- 
imuthal radio brightness profile for the chosen aspect an- 
gles (angle between the ISMF and the line of sight) are the 
isotropic and quasi-perpendicular particle injection models. 
Therefore, they conclude that the most likely orientation of 
the ambient magnetic field at the galactic latitude of SN 
1006 is parallel to the radio arcs, i.e. perpendicular to the 
Galactic Plane, in contrast to the expected direction. 

In this work, we follow a similar approach as 
IPetruk et al.l (|2009l ). but employing a 2D axisymmetric 
MHD numerical simulation and solving for the magnetic 
field rather than making a-priori assumptions on its final 
configuration. This work has a two-fold objective: to explain 
the observed morphology in radio and X-rays, and to infer 
a possible structure for the ambient magnetic field, which 
can decide which injection particle mechanism is the most 
important to explain the non-thermal emission (radio and 
X-ray) of SNRs. 



2 NUMERICAL SIMULATION 
2.1 Code description 

We carry o ut our simulations usin^ a a two dimensional Eu - 
lerian code (IPe Colle F.ll2005l : lDe Colle Raga" 20051 . 1 20061 ). 
which solves the ideal MHD equations, including radiative 
cooling. 

The equations are the following: 
| + V.(pv) = 

^ + V • (pvv + ptotI - BB) = 

1^ + V • ((e + ptot)^ - B(v • B)) = -n^HT) 

dB 

dt 



+ V • (vB - Bv) 



(1) 
(2) 
(3) 
(4) 



where I is the identity matrix, ptot — Pgas + /2 is the 
total pressure (thermal and magnetic), A(T) is the radiative 
cooling and the remaining symbols have the usual meaning. 
The equations (1-4) represent the mass, moment, energy and 
magnetic flux conservation, respectively. 

As radiative loss term A(T) we use a ta bulated coronal 
cooling function toakar no McCravl [197^ ) . The radiative 
cooling function is switched off at temperatures below lO^K. 

The constrained transport (CT) method (e.g., iTothl 

I2OOOI ) is used to conserve V ■ B — to match accuracy. 

The MHD Riemann solver uses a second-order Runge- 
Kutta method for the time integration and a spatially 



second-order reconstruction of the primitive variables at the 
interfaces (except in shocks where it reduces at first order). 
The fluxes are calculated using the HLL method. 



2.2 Initial conditions 

All calculations where carried out on a 2D axisymmetri- 
cal grid with a physical size of 12 pc and 24 pc in the 
r— and z— directions respectively, with a spatial resolu- 
tion of 1.5 X 10~^ pc. A uniform interstellar medium of 
no = 5 X 10"^cm~^ (Acero, Ballet & Decourchelle 2007) 
with a uniform magnetic field of 2/jG parallel to the z axis. 

The initial conditions consist of a sphere with radius 
0.65 pc, located at the position {ro,zo) = (0,12) pc, con- 
taining a total ejecta mass of = 1.4M0, which is con- 
sistent with a typical type la supernova. An inner region 
with constant density containing 4/7 of the ejecta mass was 
imposed, while the outer region, containing the rest of the 
mass, has a density pro file following a p ower- law density dis- 
tribution p (X fsee ljun et al.|[l996l ). This initial density 
distribution is adequate for type la SNe. 

To estimate the i nitial explosion energy w e employed 
the following equation (iTruelove k McKedflQQOl ): 



: 8.87xlO'^L>2.2 



5/3 



'^yr 



6.05x10" 



Rlno 



5/6 

+ 0.312 



(5) 



where E^i is the initial energy in units of 10^^ erg, 7^2.2 = 
D{kpc) /2.2 kpc is the distance to the object, Rb is the SNR 
observed radius in units of pc, and tyr is the SNR age in 
years. For the case of SN 1006, E^i = 2.05 considering 
Rb = 9.5 pc and tyr = 1000. In the simulations, 95% of 
the energy was set as kinetic energy and the remaining as 
thermal energy. 

Three models were used for the simulations. In model 
Ml, the remnant evolves in a uniform ISM with number den- 
sity no. In the other models, and in order to simulate the 
thermal X-ray emission as well as the deformation of the NW 
region of SN 1006, we included an ISM-density discontinuity 
at 8.2 pc from the SNR-centre, with a density three times 
higher than the ISM density. Two different directions for the 
interstellar magnetic fields were considered by changing the 
position of the density discontinuity as shown in figure [T] 
In model M2 (M3) the density discontinuity is orthogonal 
(parallel) to the symmetry axis (we recall that the symme- 
try here is referred to the numerical simulation and not the 
bright radio and X-rays arcs of SN 1006). Due to the ax- 
isymmetry, the actual density discontinuity configuration is 
a "cylindrical cloud" around the SNR. 



2.3 Synchrotron emission simulation 

From the numerical results, we calculated the synthetic radio 
emission maps by integratin g, along the line of sight, the 
radio emissivity given by fsee lGinzburg Svrovatskiill 19651 ) 



(6) 



with K being the normalisation of the electron distribution, 
B± the magnetic field component perpendicular to the line 
of sight, a is the synchrotron spectral index (which is related 
with the index of the electron energy distribution N{E) oc 
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Figure 1. Scheme of the initial numerical setup. The left panel shows the setup for model M3, where the magnetic field is parallel to 
the density discontinuity, and the right panel displays the setup for M2 where the magnetic field is perpendicular to the discontinuity. 
Notice that due to the axisymmetry assumption, the magnetic field has to be parallel to the symmetry axis and, for model M3, the 
density discontinuity is actually a cylinder about the SNR. Model Ml does not include the density discontinuity. 

lO' lO' 10^° 10^^ 10^' 

TTn| — I I — I I iiiiii | — I I m i n i 





-10 -5 



10 



Figure 2. Synthetic synchrotron maps for models Ml to M3 (from left to right) and with quasi-parallel (top row), quasi-perpendicular 
(middle row), and isotropic (bottom row) relativistic electron density distributions. The logarithmic colour scale represents the syn- 
chrotron fiux in arbitrary units (same scale is used for all models), and the horizontal and vertical length scales are given in units of pc. 
The direction of the ambient magnetic field is always along the symmetry axis towards -\-x. The vertical bright lines observed in M3 are 
artifacts introduced by the process for obtaining the simulated radio maps from a 2D axisymmetrical simulation. 



E'^, by a = (C - l)/2) whose value is 0.6 for SN 1006, and 
v is the radio frequency. 

The bilateral morphology observed in SNRs, such as 
SN 1006, can be generated by the fact that the efficiency 
of the particle acceleration at SNR shock front has a sys- 
tematic dependence on the obliquity angle between 
the normal shock and the ISMF. The literature bifurcate in 
two main lines that argue for and against the quasi-parallel 



{(t>Bn = 0°) or the quasi-perpendicular ((/)sn = 90°) mod- 
els as the most efficient mechanism for accelerating particles 
at the remnant shock front. The quasi-parallel mechanism 
i s associated with the class ical diffusive shock acceleration 
(jBlandford Eichlerl [ToSTh . while the quasi-perpendicular 
one (also called "shock drift mechanism") can be a faster 
process because of particle accelerati on by electr ic fields 
along the shock front (see for example I JokipiilflQSTl ) . 
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To include all most representative possible cases, in 
this work we study: the isotropic (no (/)sn depen- 
dence), the quasi-parallel (0sn = 0°), and the quasi- 
perpen dicular ((l)Bn — 90°) p article injection models. Fol- 
lowin g lOrlando et aD (|2007l ) (also iFulbright k Reynolds! 
I199Q| ). the isotropic case was obtained with Eq. [6] and K 
being a constant, whereas for the quasi-parallel and quasi- 
perpendicular cases, the factor K is proportional to cos^(j)Bs 
and sm^^siB? respectively. In order to identify the obtained 
simulated maps, we have added to Ml, M2, and M3, the 
words "iso", "par", and "per", for the case of isotropic, 
quasi-parallel and quasi-perpendicular particle injection. 

2.4 Thermal X-ray emission 

In order to compare with observations, the X-ray emission 
coefficient ju{n^T) was integrated along the line of sight to 
generate synthetic emission maps. 

The emission coefficient was calculated as ju{n,T) = 
in?ijj{T) (considering the low-density regime), where n and 
T are the electronic density and temperature distribution 
obtained from our numerical simulations, while '0(T) is a 
function that changes smoothly with T. '0(T) was computed 
over the energy range [0.3-5] keV using the chianti atomic 
database El and its associat ed idl software (Dere et al.l 
Il997l : iLandi k PhillipsI l20Q6l). The interstellar abs o rption 
was added employing the [Morrison McCamnionI (| 19831 ) 
model, and ta king into account a column density o f Nh = 
7 X 10^°cm"^ (|Acero et al.ll2007l : IPubner et al.ll2002l ). It was 
also assum ed that the gas is in coronal ionization equilib- 
rium flEQ. lMazzotta et al.|[T998l fl. 



3 RESULTS & DISCUSSION 

Figure [2] displays a comparison of the synthetic radio con- 
tinuum maps among all models, where in Ml the SNR is 
expanding in a uniform ISM (left column), in M2 the SNR 
collides with a density discontinuity (at a distance of 8.2 pc 
from the centre of the SNR) perpendicular to the symmetry 
axis of the simulation emulating a magnetic field perpendic- 
ular to the Galactic Plane, and in M3 the discontinuity is 
set at the same radial distance (8.2 pc) but parallel to the 
axis of symmetry, emulating a magnetic field parallel to the 
Galactic Plane (the higher density flat cloud, shown in Fig. 
1, is always parallel to the Galactic Plane). The collision 
between the SNR blast wave and the density discontinuity 
occurs at t :^ 750 yr (in the simul ation), which is in agree- 
ment with previous estimates fe.g. lAcero et al.ll2007l V 

By comparing the simulations with observations, we 
were able to discriminate some models, keeping M3par, 
M2per, and M2iso (figurefS]) as possible ones. Models M2par, 
M3iso, and M3per were left aside because the flattened 

^ (pBs the angle between the normal shock and the post-shock 
magnetic field 

^ The CHIANTI database and its associated idl procedures are 
freely available at: http://wwwsolar.nrl.navy.mil/chianti.html, 
http: / / www. arcetri . astro, it / science / chianti / chianti . html, 
htt p : / / www . damt p . cam . ac . uk/ user s / astro / chiant i / chiant i . ht ml 
^ The lEQ is an approximation for the case of young SNRs such 
as the SN 1006 giving a lower limit for the X-ray flux 



shape is parallel to the bright arcs, which is not observed 
in the SN 1006 case. The flattened shape at the NW and 
the formation of the two bright arcs are observed in M3par, 
M2per, and M2iso cases. Models M2per and M2iso show 
a very similar brightness distribution as well as intensities, 
while model M3par shows an almost two orders of magni- 
tude lower intensity for the maximum at the arcs in com- 
parison to the maximum intensities obtained for the quasi- 
perpendicular and isotropic injection models. According to 
the position of SN 1006 in the Galaxy, one would expect 
a magnetic field perpendicular to the arcs (parallel to the 
Galactic plane) , namely, our models M2per and M2iso could 
be e xcluded from the analysis if this was the case. 

iRothenflug et al.1 (|2004l ) employ a criterion which 
is basically the ratio between the flux in the interior half 
of the SNR and the flux in the limbs. They obtained an 
Rtt/s 0.7 for the radio-continuum emission of SN 1006, 
and ^ 0.3 and ^ 0.1 for the non-thermal X-ray emission 
of the remnant in the [0.8 - 2.] keV and [2.0 - 4.5] keV 
bands, respectively. IRothenflug et all (|2004h conclude that 
the bright regions observed at the limbs are probably po- 
lar caps and not an equatorial belt. They also concluded 
that acceleration of relativistic electrons proceed preferen- 
tially where the magnetic field is parallel to the shock speed. 
We computed the ratio Rj^/s for synthetic radio-continuum 
maps for the three chosen models (see figurefS]), and obtained 
0.7, 0.72, and 0.14 for M2per, M2iso, and M3par, re- 
spectively. These maps were generated considering an aspect 
angle of 90° (i.e. the ambient magnetic field lies on the plane 
of the sky). The M3par result exhibit a high contrast in ra- 
dio emission between the central region and the polar cap 
regions, which does not agree with observations. This high 
contrast increases up to 3 if an aspect angle of 11° (i.e., the 
symmetry axis is tilted by 79° with respect to the plane of 
the sky) is considered when generatin g the simulated map . 
Such aspect angle was calculated bv IPetruk et al.l (|2009h . 
The values of Rj^/s for models M2per and M2iso remain 
unchanged if this criterion is applied to synthetic maps of 
these models performed with the aspect angles obtained by 
Petruk et al. (2009) (64° and 70° for M2per and M2iso, re- 
spectively). In other words, the Rtt/s criterion is quite robust 
for both models. 

Like IPetruk et af] (|2009l \ our model is based on the 
classical MHD Rankine-Hugoniot jump conditions for the 
parallel and perpendicular components (with respect to the 
normal shock) of the ISMF. This could not be valid if the 
preshock magnetic field was randomized by an efficient non- 
linear particle diffuse acceleration. In this case the u nper- 
turbe d ISMF ca n be amplified by or ders of magnitude (jBelll 
l2004h . However. [Petruk et al.1 (|2009l ) found that, in spite of 
this, their predicted aspect angles remain almost unchanged 
for quasi-perpendicular and quasi-parallel injection models, 
after considering that the compression of a turb ulent ISMF 
does not depend on the obliquity (Bohm limit IVolk et aD 
l2003l ) and that shocks at different initial obliquity under 
magnetic field turbulent amplification become p erpendicu- 
lar inmediatly upstream ("Rako wski et al.l l2008h . The case 
of isotropic injection is discarded because it produces syn- 
thetic ma ps with no azimuth al flux variations, which is not 
observed (IPetruk et al.ll2009l ). 

The simulated thermal X-ray map for model M2 is dis- 
played in fig.Hl The symmetry axis (+x) was clockwise tilted 
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Figure 3. Direct comparison between observation (at 1.5GHz, bottom right) and models M3par (top-right), M2per (top-left), and M2iso 
(bottom-left). In order to improve the comparison, simulated synchrotron maps from models M2per and M2iso were tilted 130° clockwise 
in the plane of the sky, while the corresponding to model MSpar was tilted 40°, also clockwise in the plane of the sky. 



by 130° in the plane of the sky, in order to do a direct 
comparison with X-ray observations. This image was gener- 
ated considering an aspect angle of 90° . A bright filament is 
formed at the top right region due to the collision with the 
density discontinuity. This collision does not affect the mor- 
phology of the radio emission. Ray leigh- Taylor features are 
observed in the internal region, close to the contact disconti- 
nuity. Due to the axisymmetry^he RT fingers become bright 
lines that cross the SNR face Q In a 3D numerical simula- 



^ If an aspect angle different of 90° is employed, these bright lines 
and the bright filament would show an elliptical morphology 



tion, the RT features would give a clumpy appearance to the 
eject a emission close to the contact discontinuity interface 
(the contact discon tinuity densit y distribution was shown in 
a 3D simulation bv lMiceU et al.l l2Q09). 

Figure [5] displays the emission measure {EM) distribu- 
tion obtained by solving for the integral J rie uh dl, where 
dl lies along the line of sight. This figure was also clock- 
wise tilted by 130° in the plane of the sky and the aspect 
angle is 90°. The simulated EM can be compared to the 
parameter norm defined as 10 EM Q/4tt, where Q is 
the solid angle of the emitting gas (|Acero et al ][2003). We 
obtained an average value for the bright X-ray filament re- 
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Figure 4. Simulated thermal X-ray emission corresponding to 
model M2. The logarithmic color-bar represents the X-ray flux in 
units of erg s"-*^ cm~^ ster"-*^. The horizontal and vertical axes 
are given in units of pc. 



10 



-5 



-10 




10 



18 



- 10 



17 



-10 -5 5 10 

Figure 5. Simulated EM maps corresponding to model M2. The 
logarithmic color-bar represents the EM in units of cm~^. The 
horizontal and vertical axes are given in units of pc. 



fflon of 5 X 1 ^cm ^ which doubles the value reported by 



in IPetruk et alJ (|2009l ) . These results are reinforced when 
applying the i^^r/a criterion, since models M2per and M2iso 
yield a ratio in accordance to observations (of the order 
of 0.7) even after correcting for the aspect angle given by 
IPetruk et al.l 1^009). At the same time, the model that could 
explain the morphology observed with a magnetic field par- 
allel to the Galactic Plane (M3par) can be rejected with this 
test for synchrotron radiation. 

Furthermore, our simulated X-ray maps show a bright 
filament such as that observed in the NW of SN1006 
(jAcero et al.l l2007l ) , and the internal emission (from the 
ejected gas close to the contact discontinuity interface) re- 
vealed the RT features. Finally, our estimates of the EM are 
i n goo d agreement with the values reported by I Acer o et al.l 
confirming that the collision scenario, while capable 
of reproducing the observed thermal X-ray emission, does 
not affect the synchrotron emission. 

In summary, our results show that around of the SN 
1006 remnant, the ISMF direction is mainly perpendicu- 
lar to the Galactic plane, in agreement with previous work 
(|Petruk et al.l lioool ). Based in this fact, the best fit with 
radio observations of this remnant is achieved if the injec- 
tion of the relativistic particles at the forward SNR shock 
wave is isotropic or quasiperpendicular. This result could 
be different if nonlinear diffuse acceleration process takes 
place at the SNR shock front. In this case the model based 
on isotropic particle injection is no longer valid because the 
predicted synthetic radio maps would show a complete shell 
morphology for any position angle. Notwithstanding, models 
based on quasiperpendicular particle injection remain valid, 
with no changes in the predicted aspect angle (|Petruk et al I 
l2009l ). if this process is independent of obliquity angle (f)Bn 
(Bohm limit). 
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4 CONCLUSIONS 

In this work we present MHD axisymmetric numerical simu- 
lation of the remnant of SN 1006, taking into account differ- 
ent directions for the ISMF. Based on the numerical results, 
synthetic radio and thermal X-ray emission maps were con- 
structed. 

According to our results, and considering the morphol- 
ogy of the radio and X-rays synchrotron emission, the best 
fit is given by models M2per and M2iso, where the mag- 
netic field direction is perpendicular to the arcs, and there- 
fore perpendicular to the galactic plane. This result is in 
agreement with those obtained through analytical models 
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